Model for Assembly and Gelation of Four- Armed DNA Dendrimers 
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We introduce and numerically study a model designed to mimic the bulk behavior of a system 
composed of single-stranded DNA dendrimers. Complementarity of the base sequences of different 
strands results in the formation of strong cooperative intermolecular links. We find that in an ex- 
tremely narrow temperature range the system forms a large-scale, low-density disordered network 
via a thermo-reversible gel transition. By controlling the strand length, the gel transition temper- 
ature can be made arbitrarily close to the percolation transition, in contrast with recent model 
systems of physical gelation. This study helps the understanding of self-assembly in this class of 
new biomaterials and provides an excellent bridge between physical and chemical gels. 
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The "bottom-up" construction of new materials is one 
of the central aims of nanotechnology Q]. The recent 
synthesis of nanoparticle building blocks functionalized 
with specifically designed oligonucleotides has opened 
new possibilities for the assembly of networked mate- 
rials @> 0, 0- In this approach, a number of single 
strands of DNA are grafted on the surface of micron- 
sized particles. The interactions between the particles 
are controlled by the addition of DNA strands in solu- 
tion which have complementary sequences to the DNA 
strands grafted on the particles. The novelty of the ap- 
proach arises from the possibility of using the high sen- 
sitivity and selectivity of complementary strand recogni- 
tion to tune on and off the interparticle bonding. The 
use of DNA sequences to establish interparticle interac- 
tions provides an optimal choice for the construction of 
three-dimensional supramolecular assemblies, since DNA 
strands can self-assemble into long and fairly rigid helices 
based on sequence complementarity |5j. New materials 
can be designed by modulating the length of the bind- 
ing sequence, the length of the grafted single strand, or 
the number of grafted strands. DNA-decorated colloids 
potentially exhibit extremely rich behavior, since, in ad- 
dition to modifying base sequencing, the colloidal prop- 
erties also can be altered. This tremendous number of 
possible choices makes DNA-linked assemblies one of the 
most versatile and promising new soft-matter materials 
and thus calls for theoretical [(J Q and numerical studies 
of these systems. 

One of the key aims of such studies is the prediction 
of the three-dimensional self-assembled structure of these 
networked materials. It has been considered disappoint- 
ing that most DNA-decorated colloidal dispersions form 
highly disordered aggregates, instead of crystal-like struc- 
tures [HQ. These new materials are prone to form, in 
a reversible way, gels, i.e. disordered arrested states at 
low densities. In contrast to chemical gels, whose un- 
derstanding has progressed much further due to the con- 
ceptual simplicity introduced by the infinite bond life- 
time and the theoretical framework of percolation the- 



ory 0, 0| , the thermo-reversible formation of colloidal 
gels in the absence of phase separation and crystallization 
is still an open problem in soft condensed matter 0, • 
Recent studies have suggested that the generation of a 
thermo-reversible physical gel at low temperature T — 
i.e. conditions such that the bond lifetime is longer than 
the observation time and stresses can effectively propa- 
gate through the sample — is facilitated by limi ting the 
preferred number of bonding neighbors (valency) |l3lll4j . 
Such a constraint decreases the surface tension of a clus- 
ter aggregate, thereby destabilizing the phase separation 
process. Novel biomaterials in which a specific (and 
small) number of complementary DNA strands are at- 
tached to a common center |a, llfj should naturally lead 
to a limited valency, and hence are an interesting group 
of potential gel-forming systems. The possibility of con- 
trolling both the number of arms (the valency) and the 
bonding energy (via the number of complementary se- 
quences) makes these materials among the best candi- 
dates for checking the propensity to form gels and test 
recent explanations of physical gel formation 0] , as well 
as for exploring the types and properties of self-assembled 
biomaterials. 

In this Letter we introduce a model designed to mimic 
the tetrameric DNA complexes recently synthesized and 
discussed in Ref. 0, but whose bulk behavior has not 
yet been explored. The model is not designed to be chem- 
ically accurate, but should qualitatively reproduce the 
physical phenomena of DNA assembly. Using the model, 
we study the bulk behavior of a system of many tetramers 
and find that in an extremely narrow T range, the system 
forms a low-density disordered network via a thermo- 
reversible gel transition. In contrast with previo usly 
studied cases of thermoreversible gelation [12lll4Lllrlll7l |. 
the gel and percolation transitions of our model can be 
made arbitrarily close by exploiting the entropic contri- 
bution to bond formation made possible by the comple- 
mentary DNA binding. 

Each molecular unit of the model is composed of a 
tctrahcdral hub tethering four identical DNA-like strands 




FIG. 1: Snapshot of the simulation at low T where most 
base pairs are bonded. The red spheres indicate the tetrahe- 
dral core of each DNA tetramer. The blue tubes indicate the 
bonds along a single strand of the DNA, and the small green 
spheres represent attractive base-pair force sites. Hence the 
regular pairing of these green spheres shows the proclivity for 
complementary base pairs to attract each other. 

composed of eight connected monomers; we will refer to 
this molecular unit as a tetramer. The ordering of the 
bases beginning from the tetrahedral core is A-C-G-T-A- 
C-G-T; A, C, G, and T are the standard abbreviations 
for the bases of DNA. Bases of type A bond only with 
type T, and type C bonds only with type G. We choose 
this sequencing since it offers the possibility of forming 
bonds between different tetramers in which all eight sites 
along a strand are paired. All pairs of sites interact via 
a purely repulsive potential obtained by truncating and 
shifting the Lennard- Jones (LJ) potential 

V al (r) = Vu(r)-V LJ (r c )-(r-r c )^^ . (1) 

The cutoff r c = 2 1 / 6 er, where a is LJ length parame- 
ter. Neighboring monomers (those along a strand and in 
the central tetrahedral hub) are connected via a FENE 
anharmonic spring potential Vfene = — fc(i?o/2) ln(l — 
(r/i? ) 2 ), where the bond strength k = 30, and the max- 
imum bond extension Rq = 1.5, as used in ref. [l8j to 
study coarse-grained polymers. To model the character- 
istic rigidity of the DNA strands, we add a three-body 
potential of the form kg(l — cos 6), where 6 is the an- 
gle defined by three consecutive monomers. A value of 
ki = 5 allows for moderate flexibility of the strands, but 
prevents strands in the same tetramer complex from be- 
coming entangled. 

To simulate bond formation between complementary 
bases, each monomer along the strand has an additional 
"bonding" force site that carries the information about 
the base type. Attractive interactions are included only 
between the bonding sites of complementary bases. The 
bonding sites are connected to the monomer core using 
the same FENE potential that links the strands together. 
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FIG. 2: The fraction / of bonded DNA strands as a function 
of T. The figure shows / for the full system of 200 tetramers 
as well as for two isolated strands. The crossover from the 
unbonded to nearly fully bonded state is extremely sharp, 
but not discontinuous and it is well described by the two-state 
relation / = 1 - 1/(1 + e ( AE - TAS )/ T )) i w j t h AE = 3.53 and 
AS = 36.1. We confirm absence of hysteresis on reheating. 
We also performed a preliminary study on 64 tetramers which 
showed statistically identical behavior, suggesting there are no 
significant finite size effects. Hence, the crossover appears to 
be largely insensitive to system size. In order for the plots to 
be comparable, both systems have a strand density of 0.0056. 

The interactions between complementary bonding sites 
are modeled via a LJ potential as in Eq. QJ, but the 
truncation distance r c = 2.5a aa to include attractions. 
We choose a aa — 0.35er so that the bonding site is al- 
most completely contained in the repulsive shell of the 
monomer core. This choice prevents the bonding site 
from connecting to more than one complementary base. 
Interactions between non-complementary bases are also 
given by Eq. @J with r c = 2 1 / 6 cr Qa , so that interactions 
are purely repulsive. Fig. E] shows a snapshot of the sys- 
tem when many base pairs are bonded. We present our 
results in reduced units where length is in units of a, time 
in units of 0\J m/e, T in units e/feg (&b is Boltzmann's 
constant), and entropy is in units of fee- 

We simulate the model via molecular dynamics calcu- 
lations to generate data for the configurations and veloci- 
ties of the constituent particles as a function of time. We 
use this information to examine thermodynamic, struc- 
tural, and dynamic properties. We study a system of 
N = 200 tetramers (a total of 13600 force sites) in a 
box of length L = 52.41, resulting in a molecular number 
density n = 1.39 x 10 -3 . For this density, the approxi- 
mate average distance between the centers of tetramers 
I = n -1 / 3 = 8.96. This separation is ideal for the forma- 
tion of networks, since the distance between two bonded 
cores is 9, including the 8 base pairs in the strand and 
the tethering monomers at the core. Each simulation is 
performed at a fixed density and T, and we control T us- 
ing the Nose-Hoover method At each T we simulate 
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FIG. 3: Mean-squared displacement (r 2 ) of the tetramer cen- 
ter of mass. The dotted line indicates the asymptotic linear 
behavior of (r 2 ) that is expected when tetramers are diffusive. 
The first 4 decades in time are not shown, since the dynamics 
at this time scale are trivial. From top to bottom, the T we 
study are T = 0.140, 0.120, 0.110, 0.105, 0.103, 0.101, 0.100, 
0.099, 0.098, 0.097, 0.096, 0.095, 0.094, 0.093, 0.092, 0.090. 

2 independent systems to improve our statistics. At the 
lowest T studied, our simulations extend to more than 
10 7 time steps in order to access equilibrium behavior. 
However, we point out that this lowest T may still ex- 
hibit modest aging effects, but these will not affect our 
overall conclusions. To accelerate the overall speed of 
the simulations, we use a 3-cycle velocity Verlet version 
of rRESPA multiple time step algorithm with the forces 
separated into rapidly varying bonded forces and more 
slowly varying non-bonded forces 0| . The time step for 
the bonded forces is 0.002. 

To explore gel formation in this system, we study sev- 
eral different T values. Fig. [3] shows the fraction / of 
bonded strands as a function of T. Since the attrac- 
tive well between complementary base pairs is quite nar- 
row, we say that a base pair is bonded if the energy be- 
tween the pair is negative; we consider two strands to be 
bonded if at least half of the complementary base pairs 
of two strands are bonded. Fig. [21 demonstrates that the 
range of the transition to a highly bonded state occurs 
over a narrow range 0.09 < T < 0.11 — only 2% of 
the energy of a single bond. A simple two-state model, 
in which the open and bonded states of the strand are 
attributed to an entropy change AS — 36 and a energy 
change AE = 3.5, accurately describes the T-dependence 
of /. The large AS value shows that, in this model, bond- 
ing of strands freezes between four and five entropy units 
per base. This entropic contribution is responsible of the 
nearly first-order nature of the crossover to the bonded 
state, a feature not observed in recently studied physical 
gel models 0, El E • the consequences of which will 
become apparent shortly. For comparison, Fig. [21 also 
shows f(T) for two isolated strands. While in the bulk 



system the competition between bonding arms prevents 
reaching the fully bonded state, the fully bonded ground 
state is easily reached at low T for the case of two isolated 
strands. 

We next quantify the "freezing" of the dynamic prop- 
erties expected to occur once the gel state has formed. 
The mean-squared displacement (r 2 ) of the tetramer cen- 
ter of mass (Fig. 03) demonstrates that at high T, where 
/ w Q, tetramers are able to diffuse with little hindrance. 
Over this narrow T range, (r 2 ) becomes strongly hin- 
dered. Fig. 0Ja) shows the diffusion constant D calcu- 
lated from the asymptotic behavior of (r 2 ) = 6Dt as a 
function of 1/T. When few bonds are present, the T de- 
pendence of D is very weak. The slowing of the dynamics 
is intimately connected to the formation of bonds in the 
system; in the same narrow T region where bonding be- 
comes significant, D dramatically decreases. The nearly 
linear behavior D at low T in Fig. Ufa) indicates a lim- 
iting Arrhenius behavior, i.e. \nD ~ 1/T. More impor- 
tantly, we find a linear relation between D and (1 — f{T)) A 
(Fig-EJb)). Since / is the fraction of bonded strands, / 
also equals the probability that a strand is bonded. Thus 
1 — / is the probability that a strand is unbonded, and so 
(1 — /) 4 is the probability that all 4 strands of a tetramer 
are unbonded. The quality of the comparison between 
D and (1 — /) 4 demonstrates that the variation of D at 
fixed density is controlled entirely by the fraction of fully 
unbound tetramers. This finding mirrors the behavior 
recently found in a limited-valence model for colloidal 
gels [12| . Since f(T) can be effectively described by the 
two-state model shown in Fig. [2] at low T 

D = A)[l + exp((A£-TAS)/:r)r 4 « D e 4As e- iAE ^ T , 

(2) 

which explains the observed low T Arrhenius behavior. 

The low T Arrhenius dependence of D means that mo- 
tion will technically only cease in the limit T — > 0, but for 
practical purposes, the dynamics in the T window studied 
have already become much slower than any reasonable 
computational time. We call the T at which D can no 
longer be estimated the gel temperature T gc \ — 0.092, 
with /(Tgei) w 0.87 and D(T ge i) ~ L1 x 10_5 - A 
much more restrictive definition of T gc i, based on a dif- 
fusion coefficient ten orders of magnitude smaller — i.e. 
£>(T ge i) « 10~ 15 — would only move T gcl « 0.080, a 
decrease of just 13%. 

As we increase the length of strands, both AE and 
AS will increase, since they are proportional to the num- 
ber of bonds between the strands. Combining this with 
Eq. (J2J shows that increasing AS* will cause D to de- 
crease extremely rapidly over an even more narrow range 
of T, leading to the formation of an arrested state in a 
nearly discontinuous manner. Thus the sharpness of the 
crossover to the gel state can be tuned simply by chang- 
ing the number of bases in the strand. 

In recent model systems of physical gels ^2 Q, 0] 
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biomimetic materials for delivery and release of host com- 
ponents. 

In summary, we have shown the tendency for specifi- 
cally sequenced DNA-dcndrimcrs to assemble into amor- 
phous gel structures. In doing so, we demonstrated the 
close connection between the fraction of bonded strands 
and the dynamics of this new class of materials. As 
a final comment, we recall that in the nanotechnology 
bottom-up approach individual components are de- 
signed to assume particular tertiary structures with the 
aim of self-assembly into quaternary structures. In this 
respect, understanding the gel propensity of the designed 
nanostructure is fundamental for generating three dimen- 
sional structures with desired properties. 

We thank the NSF for support under grant num- 
ber DMR-0427239; F.S. thanks also Miur FIRB and 
MRTN-CT-2003-504712. We thank J.F. Douglas and 
L.W. McLaughlin for valuable discussions, P. Filetici for 
bringing ref. to our attention and the B.U. Center 
for Polymer Studies for hosting us. 



FIG. 4: The relation of diffusion constant D to T and bond- 
ing fraction, (a) Arrhenius plot showing that D follows an 
Arrhenius law at low T. The line uses the 2-state fitting form 
from fig. |21 (b) Dependence of D on /; the line is the best 
fit between D and the raw data for 1 — /, without using the 
2-state model. We find a power law relation D = Do(l — /) 
over all T, where Do = 0.050. Hence D is controlled by the 
fraction of unbonded tetramers. 

geometric percolation of clusters does not correlate with 
dynamic arrest, since at percolation, clusters break and 
reform continuously. Thus the T dependence of D is 
not strongly influenced by the crossing of the percolation 
locus. This is in sharp contrast with chemical gels (in 
which bonds form irreversibly), where gelation and per- 
colation coincide. In our system, the close correlation be- 
tween / and dynamics suggests that for sufficiently long 
strands the percolation of the DNA network is connected 
with the system's dynamic arrest, even though bonding 
is reversible. We locate the percolation transition using 
standard algorithms to partition tetramers into clusters 
of size s, identify spanning clusters, and evaluate the dis- 
tribution of cluster sizes n(s). Approaching percolation, 
n s ~ s~ T , with r « —2.2, the theoretical value expected 
for random bond percolation ^(J- We find the percola- 
tion temperature T p ps 0.099. At this T, molecules are 
still able to diffuse but T p is very near to T gc i — 0.092. 
For a chemical gel, T p — T gc i so that percolation ideas 
can be used to understand gel formation. Given the close 
correspondence between the length of stands and the ac- 
tivation energy for D, DNA-linked colloidal particles with 
controlled functionality make it possible to tune the dis- 
tance between T ge i and T p . The possibility of control- 
ling in a reversible way the gel transition and the on-off 
character of the transition makes the DNA-gels optimal 



[1] Z.L. Zhang, and S.C Glotzer Nano Lett. 4, 1407 (2004) 
[2] C.A. Mirkin, R.L. Letsinger, R.C. Mucic, and J.J. 

Storhoff, Nature 382, 607 (1996). 
[3] V.T. Milam, A.L. Hiddessen, J.C. Crocker, D.J. Graves, 

D.A. Hammer, Langmuir 19, 10317 (2003). 
[4] C.-H. Kiang, Physica A 321A, 164 (2003); N.C. Harris 

and C.-H, Kiang, Phys. Rev. Lett. 95, 046101 (2005). 
[5] P.J. Paukstelis, J. Nowakowski, J.J. Birktoft, and N.C. 

Seeman, Chemistry & Biology 11, 1119 (2004); N.C. See- 
man Nature (London) 421, 427 (2003). 
[6] D.B. Lukatsky and D. Frenkel, J. Chem. Phys. 122, 

214904 (2005); ibid. Phys. Rev. Lett. 92, 068302 (2004). 
[7] A.V. Trachenko, Phys. Rev. Lett. 89, 148303 (2002). 
[8] P.L. Biancaniello, A.J. Kim, and J.C. Crocker, Phys. 

Rev. Lett. 94, 058302 (2005). 
[9] D. Stauffer and A. Aharony, Introduction to Percolation 

Theory (Taylor & Francis, London, 1994). 
[10] S. Torquato, Random Heterogeneous Materials: Mi- 

crostructure and Macroscopic Properties (Springer- 

Verlag, New York, 2001). 
[11] M. Rubinstein and A.V. Dobrynin, Curr. Opin. Colloid 

Interface Sci. 4, 83 (1999). 
[12] S.K. Kumar and J.F. Douglas, Phys. Rev. Lett. 87, 

188301 (2001). 

[13] F. Tanaka and W.H. Stockmayer, Macromolecules 27, 
3943 (1994). 

[14] E. Zaccarelli, et al. Phys. Rev. Lett. 94, 218301 (2005). 

[15] K.M. Stewart and L.W. McLaughlin, J. Am. Chem. Soc. 
126 2050 (2004). 

[16] J. Dudowicz, K.F. Freed, and J.F. Douglas, J. Chem. 
Phys. Ill, 7116 (1999b 

[17] C. De Michele et al., |cond-mat/0510787| 

[18] G.S. Grest and K. Kremer, Phys. Rev. A 33, 3628 (1986); 
C. Bennemann, W. Paul, K. Binder, and B. Dvinweg, 
Phys. Rev. E 57, 843 (1998). 

[19] D. Frenkel and B. Smit, Understanding Molecular Simu- 
lation. (Academic Press, San Diego, 2002). 



